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Abstract 

We develop a version of the vacancy mediated tracer diffusion model, which follows 
the properties of the physical system of In atoms diffusing within the top layer of 
Cu(OOl) terraces. This model differs from the classical tracer diffusion problem in 
that (i) the lattice is finite, (ii) the boundary is a trap for the vacancy, and (iii) 
the diffusion rate of the vacancy is different, in our case strongly enhanced, in the 
neighborhood of the tracer atom. A simple continuum solution is formulated for this 
problem, which together with the numerical solution of the discrete model compares 
well with our experimental results. 
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1 Introduction 



Diffusion is one of the most commonly observed stochastic processes. Random 
walks, the paths of diffusing particles, are among the first typical applica- 
tions in probability theory textbooks. However, the diffusing particles, or their 
paths, are often invisible, and their effect can be observed only indirectly. This 
is the case for example in vacancy mediated tracer diffusion, when on a lattice 
of atoms a diffusing vacancy displaces atoms along its path. The process is 
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observed by following a labeled (or tracer) atom, whose steps are slaved to 
the motion of the vacancy. 

In this paper we treat theoretically the physical system of diffusing In atoms in 
a Cu(OOl) surface layer, described in detail in paper I [1]. In section 2 we review 
the hterature about tracer diffusion. We describe our discrete model in section 
3, the classical tracer diffusion problem adapted to better suit our experimental 
system. Finally section 4 gives a simple continuum formulation whose solution 
allows one to fit the experiments without adjustable parameters. 

We start with summarizing the experimental observations [1-3]. On a Cu(OOl) 
surface with up to several hundred atomic spacing wide terraces, a few surface 
Cu atoms were substituted with In atoms. Then the area was periodically 
imaged with a scanning tunneling microscope (STM), and the position of the 
In atoms was followed in time. 

To our initial surprise, the In atoms do not move in a typical diffusive way, 
but instead are stationary for some time, and then suddenly (unresolvable 
by STM) they jump to a nearby lattice site, often up to 5 or more nearest 
neighbor spacings away from their original position. These long jumps tend 
to happen at the same instant for different In atoms in the imaged area. 
To explain this, we found that the only mechanism permissible on physical 
grounds [1] is vacancy mediated tracer diffusion. 

Based on Embedded Atom Model calculations [4,5] we expect [1] that surface 
vacancies at room temperature have a low concentration, of the order 10~^ 
on Cu(OOl). They both are created and recombine at steps, and stay in the 
top layer of a terrace, as it is energetically unfavorable to dive deeper. Prom 
the typical jump rate of 10^ Hz and terrace widths of a few hundred atomic 
spacings, their lifetime is estimated to be at most of the order of milliseconds. 
(Vacancies reaching the middle of the terrace have the longest expected life- 
time.) As the long jumps of In atoms are the effect of a single vacancy, this 
short lifetime explains why the dynamics of the jumps can not be resolved with 
STM (with imaging rate up to 10 Hz), while the long waiting time between 
jumps, of the order of 10 s, enables one to distinguish independent vacancies. 

The measured waiting time between the jump events has an exponential dis- 
tribution, supporting the explanation of the jumps as independent events. The 
distribution of the jump vectors, which has been measured, will be compared 
to the numerical model in this paper. 
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2 Tracer diffusion 

The problem of vacancy mediated tracer diffusion was considered for a long 
time [6-10]. It has been solved first for the simplest case [7], when the diffu- 
sion of the vacancy is unbiased (all diffusion barriers are equal; the tracer atom 
behaves in the same way as the other atoms), the lattice is two-dimensional 
and periodic or infinite. There is a single vacancy present, it takes a nearest 
neighbor step in a random direction at regular time intervals, and has an infi- 
nite lifetime, as there are no traps. The solution is constructed by separating 
the motion of the tracer and the vacancy. The correlation of the steps of the 
tracer atom is calculated from the probability that the vacancy returns to the 
tracer from a direction which is equal, perpendicular or opposite to its previ- 
ous departure. The distribution of the tracer atom spreads with time, and on 
an infinite lattice for large times it approaches the following functional form: 



where Kq is the modified Bessel function[ll], time t measures the number 

of steps of the vacancy, and at t = the vacancy is near the tracer atom. 
The non-Gaussian shape of the spatial distribution function is typical for 
vacancy mediated diffusion. The average displacement of the tracer particle 
diverges with time. This is a direct consequence of the fact that two is the 
marginal dimension for the return probabilities in the random walk problem. 
For higher dimensional problems, the probability that the vacancy returns to 
the tracer particle is less than unity and the average displacement of the tracer 
particle remains finite. However, for dimensions equal to or smaller than two, 
the vacancy always returns to the tracer particle and as a consequence its 
displacement diverges. 

The same problem has been solved in an alternative way for all dimensions 
[8]. From this solution one can calculate the number of tracer- vacancy ex- 
changes up to time t: in two dimensions its distribution is geometric, with 
mean (logt)/7r. The continuum version of this problem has been considered 
as well in the form of an infinite-order perturbation theory [9]; the solution 
matches the asymptotic form of the lattice model. 

In a very recent study the lattice calculations were generalized to biased diffu- 
sion [10]. The difference between the tracer atom and the substrate atoms was 
taken into account by having different vacancy-tracer and vacancy-substrate 
exchange probabilities, while the rate of vacancy steps was kept constant. 
Repulsive interaction reduces, moderately attractive interaction increases the 
spreading of the tracer distribution. 




(1) 
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Although these exact solutions are closely related to our In/Cu(001) system, 
the differences, e.g. in boundary conditions and vacancy lifetime make a direct 
comparison with experiments impossible. For this purpose we develop a model 
of tracer difTusion which includes the essential properties of the experimental 
system. 



3 Discrete model 

In this section we describe a discrete model for In/Cu(001) diffusion, solve it 
numerically, and present the results. 

Our model is defined on a two-dimensional simple square lattice of size / x 
centered around the origin. This corresponds to the top layer of a terrace of 
the Cu surface, with borders representing steps. All sites but two are occupied 
by substrate atoms. At zero time the two remaining sites are the impurity (or 
tracer) indium atom, which we release at the origin, and a vacancy at position 
(1,0). This corresponds to the situation immediately after the indium atom 
has changed places with the vacancy, e.g. for the first time. 

The only allowed motion is the exchange of the vacancy with one of the neigh- 
boring atoms. The exchange rate depends on the local environment, i.e. on 
the relative position of the vacancy and the impurity atom. This takes into 
account the effect of the lattice stress induced by the tracer atom on the 
energy landscape observed by the vacancy. Each rate was assumed to be sim- 

AE 

ply proporional to the Boltzmann factor e '"s^, where AE is the activation 
energy for the considered diffusion step and /c^T is the thermal energy at 
temperature T. The diffusion barriers of the vacancy were calculated with the 
EAM method [4]. We used the interatomic potentials proposed by Finnis and 
Sinclair [5] , where the metallic cohesion is taken into account via the second 
moment approximation to the tight-binding model. The calculation was done 
on a 15 X 15 X 15 slab of copper. The system with different relative posi- 
tions of an In atom and a surface vacancy was fully relaxed until a threshold 
force of 1.875 • 10~^^ N on the atoms. The barriers were taken at the point 
halfway between the stable sites. See figure 1 for the values for the barriers. 
Since the barrier differences are large compared to /cbT, the difference in jump 
probabilities are extremely large (see below for typical values). 

When the vacancy reaches the perimeter of the lattice, its random walk is 
terminated, corresponding to the physical process of its recombination at sur- 
face steps. During its lifetime, the vacancy displaces atoms along its path, 
many of them multiple times. Thus also the tracer atom can end up displaced 
from its original position at the time of recombination of the vacancy. Av- 
eraged over the random walks of many independent vacancies, this yields a 
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Fig. 1. Special vacancy diffusion barriers near an In atom, calculated with the EAM 
method. The arrows denote the motion of the vacancy. The following values for the 
barriers were used in the calculation: Ei = 243 meV, E2 = 671 meV, Es = 503 
meV, E4 = 529 meV, E5 = 589 meV, Eq = 382 meV, E7 = 544 meV, Es = 549 
mcV, Eg = 577 meV, £"10 = 534 meV, En = 589 meV, E^ = 576 meV. For all 
other barriers E'far = 588 meV. 

probability distribution of the different displacement vectors that the tracer 
atom can make as a result of its encounter with a single vacancy. Due to the 
boundary conditions, introduced by the finite size of the lattice and due to 
the distribution of exchange rates, an analytic solution to this problem is no 
longer possible. 

When treating the above model numerically, we separate the motion of the 
vacancy and the tracer atom, as has been performed also in some of the ana- 
lytical treatments referred to in section 2. In our case of a finite lattice, this 
separation introduces an approximation, which is valid only if the tracer atom 
is relatively close to the middle of the lattice. First, we calculate the probabil- 
ities that the vacancy, released at one atomic spacing from the tracer, returns 
the first time to the tracer from equal {Peq), perpendicular (pperp) or opposite 
(Popp) directions; we also calculate the probability of its recombination (prec) 
at the perimeter instead of returning to the tracer. Knowing these return and 
recombination probabilities, we turn to the motion of the tracer atom, which 
performs a biased random walk of finite length. The direction of each step 
with respect to the previous one, and the probability that this was the last 
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step, are obtained from the return and the recombination probabihties. With 
this method we lose all temporal information about the random walk of the 
tracer, but the individual steps of the tracer atom are orders of magnitude 
too fast to be resolved with the STM, and therefore the STM observations are 
insensitive to this information anyway. 

In practice, both the return probabilities and the motion of the tracer atom are 
obtained with direct evaluation of probabilities, which has better convergence 
properties than Monte-Carlo-type methods. 

As an illustration of this enumeration method, let us consider the computation 
of the vacancy return probabilities. We assign a variable to each lattice site, 
which measures the probability that the vacancy after s atomic steps is at that 
site, while it has not exchanged with the tracer yet. Initially all probabilities 
are zero except at (1,0) where it is unity: we release the vacancy from here. 
The boundary acts as a trap for the vacancy, as well as the site (representing 
the tracer) at the origin. For each atomic step of the vacancy we update the 
site variables parallel by distributing their probability to the four neighbors 
according to the respective exchange probabilities. As the probability flows 
into the trap at the origin, we record the cumulative flow in each of the four 
directions leading to that site, which gives the return probabilities at the end. 
This iteration converges fast, and the convergence can be measured by the sum 
of the probabilities still on the lattice. The other computation, the motion of 
the tracer atom, is similar but slightly more complex. In that computation 
we assign a variable to each incoming edge of each site, which measure the 
probability that after s steps — each corresponding to a vacancy return — 
the tracer is at the given site and that it arrived from the given direction. 
In addition, each site has a variable which accumulates the probability that 
the tracer become immobile at that site. These probabilities for the tracer 
arrival and immobilization are updated iteratively according to the previously 
obtained vacancy return probabilities. 

We first tested our model on the case of unbiased vacancy diffusion, which 
would correspond to infinite temperature. The shape of the tracer distribution 
was similar to the experimentally observed one, but to achieve quantitative 
agreement the only remaining parameter of the model — the lattice size I — 
had to be tuned to astronomical sizes. This clearly shows that taking equal 
barriers is an unrealistic oversimplification for the In/Cu(001) system. 

Using the EAM barriers, for typical parameters T — 320 K and I — 401 the 
values for the return probabilities were Peq = 1 — 2.4 x 10~^, Pperp = 1-1 x 10^^, 
Popp = 4.2 X 10~^, and the recombination probability prcc = 1-1 x 10~^. These 
values depend weakly on / [e.g. the dependence of the mean square displace- 
ment, calculated later in this paper from the return probabilities, is logarith- 
mic: (r^) oc log(///o)]- This is a consequence of the fact that two is the marginal 
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dimension for the return problem of the random walker. In higher-dimensional 
space the vacancy does not necessarily return, the return probabilities are 
asymptotically independent of the lattice size, and the final distribution of 
the tracer — also independent of lattice size for large lattices — takes the 
following form, e.g. in three dimensions [9]: p(r) = Ci5(r) + C2r~-^ exp(— Ar). 

For the case of an indium impurity in a Cu(OOl) lattice, the diffusion barrier 
for a vacancy exchange with the indium atom is considerably lower than all 
other barriers. Therefore, in most cases the vacancy returns from the direction 
of its previous departure, and the individual moves of the tracer atom are 
strongly anti-correlated. Both the numerical and the theoretical treatment 
are simplified significantly if we do not have to follow the large number of 
ineffective "back and forth" exchanges of the vacancy and the tracer atom. 
For this purpose, consider the small probability e that the vacancy does not 
return from the same, equal direction. Thus Peq = 1 — e. In the case of indium 
in Cu(OOl) the EAM calculations yield e — 2A ■ 10"''' (see above). We now 
define 



Pperp Pperp/^ 
Popp Popp/^ 
Prec Prec/^ 



and have 



2pperp ~l~ Popp ~l~ Prec 



(2) 



(3) 



If we represent the quasi-bound state of the rapidly exchanging (on average 1/e 
times) vacancy and impurity atom with the position of a bond of the original 
lattice, then the vacancy-tracer pair walks on the bonds of the original lattice. 
The pair steps to each of the four perpendicular bonds with probability Pperp/2, 
and to each of the two parallel ones with Popp/2, see figure 2. (The factors 1/2 
reflect the probability for the vacancy to escape either at the right or the 
left side of the quasi-bound position, which is | for both sides in the limit of 
vanishing e. The advantage of this approach is twofold: the path of the tracer 
is made of fewer steps (beneficial for numerics), and the bond-to-bond steps 
are now independent (beneficial for theoretical treatment). 

Using this, the tracer atom described as if it forms a pair with the vacancy on 
one of the bonds adjacent to its original site, walks on the bonds lattice, and 
at the end (which happens with probability p^ec after each step) it is released 
at either end of the last visited bond. Results for the probabilities of the 
different jump lengths (beginning-to-end vectors of these paths) are shown on 
figure 3. Note, that the model calculations in figure 3 contained no adjustable 
parameters. 
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Fig. 2. a) The bond-to-bond steps of the vacancy-tracer pair. If we view the va- 
cancy-tracer pair to Uve on the center of each bond, indicated by the labels, then 
the pair hops between these sites. The atomic lattice is drawn with thin lines, b) The 
lattice (made of bonds of the atomic lattice) on which the vacancy-tracer pair walks. 
The neighbors that can be reached by a single step from the site at the center are 
shown with thick lines. This lattice can be considered as a square lattice, rotated by 
45° with respect to the atomic lattice, with some extra diagonal bonds. Note that 
this is not a Br avals lattice. 




I I I I I I I I I I I I I I I I 

1 2 3 4 5 6 7 

jump length r [nearest neighbor spacing] 

Fig. 3. The probabilities of the jump lengths of the tracer atom for T = 320 K and 
I = 401 lattice spacings. Filled circles correspond to experimental values (measured 
at this temperature and terrace size), open circles are from the model described 
in the text, and the solid curve is the continuum solution described in section 4. 
The data shows no significant directional structure: the dependence on the length 
of the jumps is monotonic with good approximation. (Each dataset is normalized 
separately such that the probabilities corresponding to a subset of the jump vectors, 
1 < |r| < 6, add up to unity. These are the probabilities that are determined with 
good accuracy in the experiment.) 
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Fig. 4. The temperature dependence of the probabihty that a tracer atom had an 
encounter with a vacancy but its net displacement was zero, plotted for two different 
lattice sizes. The dependence is weak. 

A general advantage of numerical modeling is that we have access to quantities 
which are difficult or impossible to measure experimentally. One example in 
our simulation is the probability that a tracer atom had an encounter with a 
vacancy, but its net displacement was zero. The temperature dependence of 
this quantity is plotted on figure 4. Since the dependence is weak, the assump- 
tion in the experimental measurements to associate In-vacancy encounters 
with visible (non-zero) jumps of the In atom is justified. 

For a given set of diffusion barriers, our model has two parameters: the tem- 
perature and the lattice size. When we compare results with experiments, both 
can be independently obtained, and in principle there are no adjustable pa- 
rameters. For example the case of T = 320 K (figure 3) — using the distance 
to the nearest step in the experiment as lattice size in the numerical model — 
gave a good match with the experiment, but for the measurements at other 
temperatures we had to adjust the lattice size to obtain good agreement. Al- 
though the best fit lattice size in some cases was a factor of 2-3 smaller than 
the measured distance to the nearest step, the change in the return proba- 
bilities was much smaller, as they depend logarithmically on the lattice size. 
Undetected defects — acting as trap for the vacancies — closer to the In atom 
than the nearest step could also be accounted for this difference. 



4 Continuum model 



In the previous section we solved numerically the version of the tracer diffusion 
problem which was relevant for the In/Cu(001) experiments. Although this 
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already enables full comparison, the numerical solution has the disadvantage 
that it cannot be described with a few parameters. In this section we develop 
a simple continuum description, where the overall shape of the jump length 
distribution is described with a single parameter. 

We use our previous results for the return and recombination probabilities of 
the vacancy, and consider the random walk of the tracer-vacancy pair on the 
bond lattice. Let ^(r, n) denote the probability that the tracer-vacancy pair is 
at position r and at instance n, where n counts the number of steps the tracer- 
vacancy pair takes. Since the subsequent steps of the pair are independent, we 
can write an effective diffusion equation for the evolution of g{r,n): 

— = Deff V g-cg. (4) 

on 



The first term on the right hand side corresponds to the steps the pair takes on 
the bond lattice, here denotes the mean square displacement per step of 
the pair. The second term corresponds to the recombination of the vacancy^; 
at this point the pair breaks up. In the continuum approximation for space 
and n, the solution for a Dirac-delta initial condition at the origin is 

^(r, n) = — ^ exp (--^ cn] . (5) 



The final distribution of the tracer atom after the vacancy recombined is 
obtained by the integration of the loss term in Eq. (4): 



p(r) = / cg{r, n) dn = — -r— Kq 



r 



(6) 



where Kg is the modified Bessel function of order 0. The functional form of 
this solution is similar to that of the infinite lattice system in section 2, in spite 
of the differences introduced by the finite vacancy lifetime and the different 
boundary conditions. 

The mean square displacement (r^) is directly proportional to the square of the 
width of the Bessel solution, apart from lattice corrections. We can determine 
the width of the Bessel solution from the parameters D^s and c, which are 

^ Note that since n counts the number of steps of the vacancy-tracer pair on the 
bond lattice, the term —eg does not imply that the vacancy can recombine at any 
lattice site — in fact it recombines at terrace steps, between subsequent returns to 
the In atom. 
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obtained from the return probabilities in the Appendix: 

This continuum solution is shown in figure 3. It closely follows the numerical 
solution of the model, even for relatively small distances from the origin, where 
one would expect stronger lattice effects. 



5 Summary 

In this paper we described a model for vacancy mediated tracer diffusion on a 
finite lattice with absorbing boundaries for the vacancy. These boundary con- 
ditions were appropriate to model the vacancy mediated diffusion of In atoms 
embedded in the top layer of Cu(OOl) terraces. In addition to the numeri- 
cal solution of the discrete model, we set up a simple continuum formulation 
of the model. The spatial distribution of the tracer atom in the continuum 
solution has a modified Bessel function profile. This form of non-Gaussian 
distribution is typical for tracer diffusion assisted by other diffusing particles. 
In order to enable a quantitative comparison with the STM measurement of 
the In/Cu(001) system, we introduced modified vacancy diffusion rates near 
the In atom, calculated with the EAM method. The modified rates affect the 
width of the Bessel function, without changing the functional form of this 
characteristic distribution. 
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A Obtaining D^q and c from the return probabilities 



In this Appendix we calculate the parameters of the continuum model, Defr and 
c, from the return and recombination probabilities in Eq. (2). The probability 
Po^"+^^ that the vacancy-tracer pair is at site "0" (not necessarily the origin) 
after n + 1 steps of the pair is contributed from sites 1, 2, . . . , 6 (see figure 2a; 
no upper index refers to probabilities after n steps): 

Pt^'^ = ^{Pi + P2 + P3 + P.) + ^-^{P, + P.) . (A.l) 



During one step in n, the change in Pq is (using Eq. 3) 
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Po("+i) - Po = ^(Pi + P2 + Ps + P4 - 4Po) 

+^(P5 + Pe - 2Po) - AecPo . (A.2) 

The finite difference on the left hand side approximates the derivative, and the 
terms on the right hand side can be collected to form a lattice approximation 
to the Laplacian: 

V'P = aZ^JPi + P2 + P3 + P4- 4Po) (A.3) 



and 



dlP = aopp(A + Pe - 2Po) ^ ^V^P , (A.4) 



where the distances Opcrp = 1/ and Oopp = 1 are in units of nearest neighbor 
spacing of the atomic lattice. In the approximation of Eq. (A.4) we assumed 
on average equal second derivatives in all directions. After these substitutions 

^^ Pperp+Popp ^2p_g^p^^ (A.5) 
^ » ' c 



and the coefficients Des and c can be read off. 
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